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Abstract 

This paper presents a study of the two-predators-two-preys discrete-time Lotka-Volterra model with self- 
inhibition terms for preys with direct applications to ecological problems. Parameters in the model are 
modified so that each of them has its own biological meaning, enabling more intuitive interpretation of 
biological conditions to the model . Moreover, the modified version is applicable to simulate a part of a 
large ecosystem, not only a closed predator-prey system with four species. An easy graphical method of 
analysis of the conditions on parameters ensuring long persistence under coexistence is presented. As an 
application, it is explained that a predator specie who feed on relatively small number of preys compared 
to the other predator species must be selective on the available preys in order for long persistence of the 
ecosystem. This may be regarded as a theoretical explanation of the existence of flush-pursuer birds, those 
who uses highly specialized hunting strategy and cross-adapts to the ecosystem relative to the ordinary bird 
species. 



1. Introduction 

This paper is a study of four-dimensional discrete-time Lotka-Volterra model by the method of lineariza- 
tion and simulations, with application to a realistic ecological system. The Lotka-Volterra models, proposed 
by Lotka [1] for the description of chemical reactions and by Volterra [2] for the description of interacting 
populations, have been the starting point in many dynamical models. The original and most commonly 
studied version of this model describes a continuous-time dynamics, whereas in real ecosystem the changes 
in populaions of each species due to competitive interaction cannot occur continuously. Hence it would be 
more suitable to use a discrete-time dynamical system, represented by difference equations or, more properly, 
by the iterated application of maps, are often a more suitable tool for modeling the dynamics in ecosystem. 

The choice of discrete-time dynamic models is crucial because they may exhibit more complicated dy- 
namic behaviors then continuous-time models. Indeed, even in one-dimensional discrete dynamical systems 
represented by iterated quadratic maps, like the well known logistic map, periodic and chaotic trajectories 
can be easily observed [3]. Similarly, two-dimensional discrete dynamical systems, in particular discrete time 
Lotka Volterra equations, can exhibit a wide variety of complicated asymptotic behaviors, from convergence 
to a fixed point or a periodic cycle until quasi-periodic motion along a closed invariant curve or even an 
erratic motion inside a bounded chaotic attractor (see e.g. [15] , [6], [7]). 

The local asymptotic behaviour near a fixed point of the system can be analyzed by linearization, which 
is basically an approximation of the system by a linear system defined by the Jacobian at the given fixed 
point. The global behaviour of a discrete-time dynamical system, however, is far more interesting and 
complicated. One of a few methods available for the study of it is the elegant theory Lyaponov; if one 
can find a scalar function that satisfies certain conditions, called the Lyaponov function, then one can say 
about the global asymtotic behaviour to certain points or bounded in certain region. But since their is no 
systematic way to find the Lyapunov function of a given dynamical system, this is often a hard task. The 
method of critical sets is another way of studying the global behaviour of discrete-time dynamical system 
which is applicable if the transition map is non-invertible, as it happens to be in the case of Lotka Volterra 
maps. Such method is a powerful tool introduced in the Seventies, by which the global dynamical properties, 
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especially the bifurcations, can be usefully characterized(see [T5], [T7])- See [T3] for such a study of global 
properties of three dimensional Lotka-Volterra map by means of critical set. But the method of critical set 
applied to a map of dimension greater than two, may contain a very hard and challenging task, due to the 
difficulties met in the computer-assited graphical visulaization(see e.g. [TH], Q35]) There is also an excelent 
mathematical analysis of a special kind of three-dimensional Lotka-Volterra map by Gardini et al. [2"D] 

However, those methods are too difficult to be practical for biologists who wants to use the mathematical 
modeling to an actual biological problem in order to support an emperical hypothesis or make a theoretical 
hypothesis to guide further emperical study. To present a sound and practical modeling stretage which 
is easily applicable to some interesting biological problems, we need to satisfy the following two points. 
First, the mathematical model must reasonably simulate the biological systems in interest. Second, the 
method of analyzing such a mathematical system need to be handy and could draw a meaningful biological 
implications. In the present paper, we concern a four dimensional discrete-time Lotka-Volterra system with 
2 preys - 2 predators 



which we consider as the one satisfying the two requirements in some acceptable extent. First, the system 
in concern can simulate biological situations : 1) a system where a specie of generalist predator and a specie 
of specialist predator which is selective at a certain prey specie competes over two prey species; 2) a system 
with a predator specie heavily dependent on two prey species and with another predator specie which is 
less dependent one the two prey species due to the presence of other source of food; 3) a closed system with 
single prey and predator species to which a new predator specie is introduced by migration or mutation. 4) 
a closed system with single prey and predator species to which a new prey specie is introduced by migration 
or mutation. Second, as a method of analyzing such system, we introduce an application of linearization 
which enables us to draw some biologically meaningful conclusions about the global behaviour of the system 
and conditions on prey and predator species for coextence in the current system. 

2. The two-dimensional Lotka-Volterra map 

3. The three-dimensional Lotka-Volterra map 

The three-dimensional Lotka-Volterra map T3 : (xi,X2,Xs) l— > (^'i, %' 2 , x' 3 ) is defined by 



4. The four-dimensional Lotka-Volterra map 

The general four-dimensional Lotka-Volterra map T4 can be defined as similarly as the three dimensional 
map T 3 . However, we concentrate on the following specialization T' A : (x\, X2, X\, X 2 ) i-> (x[, x' 2 , X[, X£) 
defined by 



Xl {t + 1) = aBi(t)[n - feisi(t) - anXtit) - a 12 X 2 (t)} 
x 2 (t + 1) = x 2 (t)[r 2 - k 2 x 2 (t) - a 2 iXi(t) - a 22 X 2 (t)} 
Xt(t + 1) = X x (t)[-pi + bnxxit) + b 12 x 2 {t)] 
X 2 (t + I) = X 2 {t)[-p 2 + b 2lXl (t) + a 22 x 2 {t)] 



(1) 






x[ = a;i(ri - kixi - biiXt - b 12 X 2 ) 
x 2 = x 2 (r 2 - k 2 x 2 - b 2 iXi - b 22 X 2 ) 
X[ = Xi(-pi + c n xi + c 12 x 2 ) 
X 2 = X 2 (-p 2 + c 2X x 1 + c 22 x 2 ) 



(3) 
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where all the parameters are positive real numbers. The discrete-time dynamical system determined by the 
repeated application of this map T' A simulates a predator-prey ecosystem consisting of two prey and predator 
species; the lower and upper case variables denote prey and predator species, repectively. The underlying 
biological assumptions for this model are as follows. 

(a) . All species have overlapping generations; the populations of all species are updated simultaneously. 

(b) . There is no direct interaction within preys or predators. 

(c) . We does not consider any age structure in the model. Any two individuals of the same species are 

equivalent (have the same traits) regardless of their age. 

(d) . Without the presence of predators, the prey populations follow logistic growth. 

(e) . Without the presence of preys, the predator populations follow exponential decay. 

(f ) . The amount of predation of prey i by predator j and reproduction of predator j due to the consumption 

of prey i at each generation are propotional to the product XiXj of their population, as such product 
represents the number of encounters during that generation. 

Analyzing the map Q, it is easy to see that the coordinate axes, planes, and spaces are invariant sets, 
i.e., if a coordinate of x S t 4 is zero, then the same coordinate of T^x) is zero. For example, a trajectory 
starting from a point such that xi(0) = remains trapped inside the 3-dimensional space x\ — and is 
determined by the restriction of ^ obtained from setting x\ = 0. One easily realizes there are total 15 
fixed points with at least one vanishing coordinate. Now suppose is a non- degenerate fixed point of ([3]), 
i.e., a fixed point without zero coordinate. It is clear that such a fixed point is given as a solution of the 
four linear equations obtained form canceling out the degree 1 factors in each of the four equations in ([3]). 
In other words, we have 
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Of course, the 4x4 matrix inverse in Q may not exist and there could be no or many non-degenerate fixed 
points of Such cases are not of our interest here, however, and as long as we are interested in modeling 
ecosystems, there is no serious problem in assuming that such inverse matrix exists since the probability 
that the inverse matrix does not exist when the parameters of ^ are arbitrarily given, is zero. Hence we 
may assume that ^ has a unique non-degenerate fixed point given by Q. We will be more interested in 
such a non-degenerate fixed point if all of its four coordinates are positive. We may call such fixed point 
positive. 

After all, we are interested in the condition for the parameters of the map T4 such that it has a unique 
positive equilibrium and the system persists near the fixed point for a relatively long time. This is the 
situation when the corresponding predator-prey system persists with all the four species coexisting. For this 
task, we shall use the standard technique of linearization, which we will briefly summarize. 

As one studies the trend of a curve near a certain point by looking at the tangent line, one can study a 
non-linear system near a fixed point by approximating to a well-understood relatively simple system, which 
is, so to say, tangent to the non-linear system at the given fixed point. Such a simple system is called linear 
system, whose transition map is given by matrix multiplication. For example, if A be a 4 x 4 matrix with 
real entries, then we have a four-dimensional linear system given by the transition map 

x' = Ax. (5) 

Thus gien an initial condition x(0), the state x(n) after n iteration of the map ^ is given by 

x(n) = A"x(0). (6) 
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If, in particular, the initial condition x(0) is an eigenvector of A with eigenvalue A, then the trajectory will 
be 

x(n) = A"x(0). (7) 

Hence x(n) converges to the origin if |A| < 1, and diverges from the origin if |A| > 1, in both cases with 
speed |A|. If one requires all the trajectories converge to the origin, it means that the matrix A n converges 
to the zero matrix. In this case we say that A is zero-convergent and the origin is a stable fixed point of 
(§■ If some trajectory tends to infinity, then the origin is said to be unstable. More loosely, if one requires 
that every trajectory is trapped in certain region, it means that there is a positive real number M such 
that every entry of the matrix A n has absolute value less than M for all n > 1. In this case we call A 
power-bounded. Evidently, a matrix A is zero-convergent only if every eigenvalue must has modulus less 
than 1, and power-bounded only if every eigenvalue has modulus no greater than 1. It is well-known that 
these conditions are also sufficient. This is the condition for long persistence of a discrete linear system. 
Hence one can determine the dynamics of a linear system by looking at the eigenvalues of the transition 
matrix. 

On the other hand, one can find such a linear system which is tangent to a given non-linear system at 
a given fixed point by finding the Jacobian matrix of the transition map at the fixed point. In our case, 
let x* = (x\ , x 2 , ^1 , X 2 ) be a fixed point of the map T' A . By using the equilibrium condition, the Jacobian 
matrix J4 of T[ at x* is given by 



1 — k\x\ — a\\x\ — a\ 2 x\ 

1 - k 2 x 2 -a 21 x* 2 -a 22 x 2 

&11A7 b 12 X{ 1 

h 2l X* 2 b 22 X^ 1 



(8) 



The derived linear system 



x' = Jjx. (9) 

is the one to which we approximate the dynamics of our system near the fixed point x*. As we have seen, 
the system Q is locally stable at the fixed point x* if the four eigenvalues of J' 4 is less than 1 in modulus, 
and unstable if at least one eigenvalue has modulus greater than 1. Moreover, one can guess that in order 
for the trajectories near the fixed point to be bounded, the Jacobian J' 4 must be at least power-bounded. In 
other words, if the system ^ shows a long persistence near a fixed point, then it is very much likely that 
the eigenvalues of the corresponding Jacobian matrix is no greater than 1 in modulus. Even though there 
are not easy algebraic manipulations to express the four eigenvalues of J4 in terms of the parameters, the 
formulae (4) and (5) enables us to determine whether the system Q with a particular set of parameters has 
a stable non-degenerate fixed point or not. 



4-1- Modeling 

Although the map T[ is handy for mathematical analysis, it is not so practical to use to simulate a 
particular type of ecosystem; it is not clear how to assign suitable values of the parameters to reflect 
important characteristics of the subject ecosystem. Look at the equations in ^ . The parameters r and k 
describes the intrinsic population dynamics of preys, and p describes the dependency of predators to preys. 
What describes the interaction between the four species are the parameters b and c, denoting the predation 
rate of preys and reproduction rate of predators, respectively. However, these two parameters are the results 
of combined effect of various factors between predators and preys. We assume the predation of preys is 
determined by the following three factors: 

(1) hunting efficiency of a particular predator specie over a particular prey specie; 

(2) adaptation of a particular prey specie to a particular predator specie; 

(3) search rate of a particular predator specie. 
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Let eji be the hunting efficiency of predator j over prey i, aji be the adaptation coefficient of prey i to 
predator j, and Si be the search rate of predator i. We assume the number of encounters between predator 
j and prey i is SiXjXi, and the number of prey j predated by predator i during one generation is given by 

SidijejixiXj (10) 

Note here that the larger the value of a^j , the more the prey i is adapted against predator j . Hence we 
may factor the predation rate b{j by bij — SiCiijeji. On the other hand, the amount (10) of hunted prey i 
will be consumed and contribute to the reproduction of predator j for the next generation. We multiply 
another coefficient qji, the conversion ratio, to (10), to measure the amount of reproduction of predator j 
due to the consumption of prey i. Hence we may factor the reproduction rate Cji by Cji = q^a^e^. Finally, 
we may write k,- L = ri/K i} where is the carrying capacity of prey i. Now we reformulate the map T' A as 
follows: 

'x[ = nxi(l - fi-) _ sianenXiXi - sia 1 - 1 ei 1 x 1 X 2 
T , . < x 2 = r 2Xi{l - j^) - s 1 a 2 ie 12 x 2 X 1 - s 2 a 22 e 22 x 2 X 2 

X[ = X 1 - p 1 X 1 + sxqnanenXxXx + s 1 q 12 a 21 e 12 x 2 X 1 
^X' 2 = X 2 — p 2 X 2 + s 2 q 21 a 12 e 21 x 1 X 2 + s 2 q 22 a 22 e 22 x 2 X 2 

This model will not only be interpreted as a model of a closed predator-prey system with four species, but 
also a model of a portion of a large ecosystem, according to the values of the parameters and their biological 
meaning. 

Figure 1 illustrates six numerical simulations that show various dynamical patterns of the system (??). 
The population of the four species are plotted versus the number of generations elapsed. Red, green, black, 
blue dots refer to Xi,x 2 ,X 1 ,X 2 , respectively. 
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Figure 1: (a), (b), (c) have common parameters r\ = r 2 = 1.5, K\ = K 2 = 10 4 , s± = s 2 = 0.01, Pi = P2 = 0.3, and 
<?n = 912 = 921 = 922 = 0.02 but different hunting rate; (a) an = ai2 = 0.3, 021 = 0.2, a22 = 0.5 (b) an = 0.2, ai2 = 0.45, 
a2i = 0.5, 022 = 0.25 (c) an = 0.24, ai2 = 0.6, 021 = 0.56, 022 = 0.27. (d) has the same value of r;s and qijS, but 
Kt = K 2 = 10 5 , pi = p 2 = 0.01, si = s 2 = 0.001, an = 0.2, a 12 = 0.4, a 2 i = 0.31, a 22 = 0.29. (e) and (f) have the same 
parameters; their set of parameters differ from that of (d) only in the a's; an = 0.4, ai2 = 021 = 0.2, 022 = 0.35. The first four 
have the same initial condition (100,100,100,100), whereas (700,2000,2300,3000) for (e) and (10,100,200,2500) for (/). The 
non-degenerate fixed points are (a) (3.3333e + 3, 1.6667e+3, 277.7778, 83.3333), (b) (1.7143e+3, 2.5714e + 3, 140.8163, 192.2449), 
(c) (1.8252e + 3,1.7699e + 3, 132.8351,162.0379), (d) (833.3333, 833.3333, 450.7576, 4.5076e + 3), (e), (f) (750.0000m 1.0000e+3, 
2.2406e+3, 2.9625c+3) 

The first four examples starts at the same initial point, and the first three cases show different dynamics 
depending on the values of hunting rates. Hence the two prey species are identical with the same intrinsic 
rate of increase and carrying capacities, and the two predator species have the same dependency on preys and 
conversion ratios. These parameters are suitably chosen to allow long persistence. However, the strategies 
of predators are what make the difference. In (a), predator 1 is equally efficient at hunting both prey 1 
and 2 with an = an = 0.3, but predator 2 is specialised in hunting prey 2, with 0,21 = 0.2 and 022 = 0.5. 
Hence it simulates a system with a generalist and a specialist predator. The system converges to the unique 
positive fixed point, that is, there is a stable coexistence at the fixed population. 

In (b), both predators are specialized in different prey species. The dynamics is a regular oscillation 
around the fixed population, and the persistence is very strong. Example (c) also simulates a similar system, 
but now the predator 1 is slightly more specialized than before. The dynamics is a irregular oscillation around 
the fixed population, with characteristic pattern that the period of oscillation continuously increases and 
decreases alternately. 

In (d) we have 10 times more carrying capacities for the preys and less dependencies of predators and 
search rates. Predator 1 is specialized in hunting prey 2, while predator 2 is almost equally efficient on 
both preys. Its dynamics is interesting, showing a domination by the generalist predator 2, with alternating 
dominance of prey 1 and 2. Under the dominance of predator 2, the number of predator 1 keeps decreasing 
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and when this number goes down below a certain level, the prey 2 suddenly bursts out and prey 1 decreases 
rapidly. Due to the abundance of prey 2 to which the predator 1 is specialized, the number of predator 1 
gradually increases. But at some point large number of predator 1 almost extinguishes the prey 2, and then 
prey 1 now thrives with rapid oscillation. 

The last two examples show the dependence on the initial conditions. Both (e) and (f ) have the same pa- 
rameters but different initial populations. The initial population for (e) is more close to the fixed population, 
and it shows almost regular oscillation around the fixed populations. On the other hand, the system of (f) 
undergoes increasing or decreasing oscillations and then suddenly collapses after about 600th generations. 

Nonetheless, all the six examples persist for a very large number of generations, mostly around the fixed 
populations. Such a long persistence around the fixed population would be not so much likely if there is 
an eigenvector with eigenvalue whose absolute value far exceeds 1, since then any component of the initial 
population vector along the eigenvector will be rapidly repelled from the fixed point. Indeed, the table below 
shows that the eigenvalues in the six examples are mostly less than 1, and some of them are slightly greater 
than 1. 

Notice in the three unstable cases, the absolute values of the four eigenvalues are close to the critical 
value 1, especially the one exceeding 1. It seems that a long persistence of the system around the unstable 
fixed point is only possible if the spectral radius of the Jacobian, i.e., the maximum value of the modulus of 
the eigenvalues, which exceeds 1, is closed to 1. The decaying oscillation toward the fixed point shown in (a) 
is due to the complex eigenvalues less than 1 in modulus. The presence of complex eigenvalues greater than 
1 in modulus as in the other three unstable cases may contribute to the persisting oscillatory behaviors. 
However, it would not be possible to predict the global pattern of the system solely from the local data of 
eigenvalues. 





Eigenvalues and their modulus 


spectral radius 


(a) 


0.8147+0.5277i 
0.9707 


0.8147-0.5277i 
0.9707 


0.6597 


0.9609 


0.9707 


(b) 


0.8389+0.6569i 
1.0119 


0.8389-0.6569i 
1.0119 


0.8397+0. 1285i 
0.8495 


0.8397-0.1 285i 
0.8495 


1.0119 


(c) 


0.8652+0.5824i 
1 .0485 


0.8652-0.5B24i 
1 .0485 


0.8652+0.19301 
0.8865 


0.8652-0.1 930i 
0.8865 


1.0485 


(d) 


0.9937+0. 1218i 
1.0012 


0.9937-0. 1218i 
1.0012 


0.993 8+0.0112i 
0.9938 


0.9938-0.0112i 
0.9938 


1.0012 


(e), (f) 


0.9933 + 0.121 7i 
1.0007 


0.9933 - 0.121 7i 
1 .0007 


0.993 8+0.0112i 
0.9942 


0.9938-0.0112i 
0.9942 


1.0007 



Figure 2: Eigenvalues and their modulus of the corresponding Jacobian matrix 



local repelling force, attracting force 

global attracting force : self-inhibition term for preys inhibit total collapse 

4-2. Stability diagram analysis 

In this section, we introduce the notion of stability diagram for the system (111. We may fix the param- 
eters other then those hunting efficiencies ay, and investigate the quadraples (an, 0,12, 0,2x1 0,22) for which 
the system (11) have stable positive equilibirum. But the result would be plotted in 4-dimensional space, 



which is impossible to visulaize. To reduce the dimensionality, we may simply assume that 



an + a i2 



1. 



a 21 + a 22 



1. 



(12) 



The biological assumption under the constraint ( 12 ) is that each predator species has a limited hunting 
ability to be distributed over two prey species. For example, if predator 1 can hunt down prey 1 at rate 0.7, 
then the hunting rate for prey 2 must be 0.3, the remainder of the former. But of course two predator species 
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can have different totai hunting ability, and this point is reflected in the search rate. Namely, if predator 
1 have search rate twice of that of predator 2, meaning that predator 1 encounters preys twice more than 
predator 2, then predator 1 actually have twice total hunting ability than that of predator 2. 



Moreover, the constraint (12 1 does not harm the generality of our model, since one can normalize the 
parameters of a given system (8) by replaing each paratmeter ay by ■ and Si by s£ , where 

, _ On , _ d2i , _ 021 , _ «22 

On — i O21 — i 21 — ' ^22 — 

an + O12 On + O12 021 + 022 O21 + 022 

s[ = Si(on +ai 2 ), s' 2 = s 2 (o 2 i + a 22 )- (13) 

Notice that the "total hunting ability" of prey i is reflected in the factor (a^i + a^ 2 ) of the normalized serch 
rate = Si(an + )• Also, now it make sense to call the hunting rate as "hunting probability" since the 
normalized hunting rate a^ has value < a™ < 1. We may call a given four-dimensional Lotka-Volterra 



map (11) normalized if it satisfies the constraint (12) 



By normalization, the tuple (an, 012, 021, 022) becomes (an, 1 — an, 1 — 022, 022), and hence determined 
by the two parameters an and 022- Thus we can plot all the pairs (an, 022) in the unit squre S = {(x, y) £ 
M 2 \ < x,y < 1} for which the given normalized system has positive and stable fixed point. We call the set of 
all such pairs (011,022) the stable region, and such diagram as the stability diagram of the given normalized 
system. As we have observed in the simulations in section 4.1, a long persistence of a given system around 
the non-degenerate fixed point might be possible only if the spectral radius of the system's Jabobian is very 
closed to 1. Even though we are dealing with a discrete dynamical system, the spectral radius of Jacobian 
for a given system is a continuous function in the parameters. In other words, if the two tuples (011,022) 
and (a'ii,a 22 ) in the stability diagram are very close, then the corrdsponding spectral radii will be closed 
as well. Hence, we obtain the following principle for long persistence: a system corresponding to the tuple 
(an, 022) would be persistant for a long periods of generations only if the tuple (011,022) is closed to the 
stable region. 



Given the values of parameters of system ( 11 ), one can readily determine whether the corresponding sys- 
tem has a positive stable fixed point algorithmically, using the formulae Q and ([8]). In stepwise discription, 

(a) Find the nondegenerate fixed point of the given system by using the formula Q and check if all the 

coordinates are positive. 

(b) Using ^ , find the Jacobian matrix of the system at the fixed point and calculate the four eigenvalues. 

(c) (an, 022) is a stable tuple if all of the the four coordinates of fixed point x* are positive, and the 
maximum of the four eigenvalues in step (b) have minimum absolute value less then 1. 

Therefore, one can easily write a program in Matlab that draws the stability diagram for given parameters 
except the hunting rates, by implementing the stepwise algorithm as described above. So this could be a 
practical and useful tool that visualizes the results of localization analysis, to "see" what combinations of 
an and 022 may lead to a long persistance of the system with respect to the other given parameters. We 
will use this method for various set of parameters that represent certain type of ecosystems observable in 
nature. 

Localist vs. Globalist. First, consider a predator-prey system where "globalist" predator and "lo- 
calist" predator compete. Globalist predators feed on various source of food since they search for preys on 
large area, whereas localist are localized so that they have very limited source of food compared to globalist. 
We are going to emulate this situation via our 4-dimcnsional model. Suppose the localist feed on two prey 
species. Then they ought have high search rate and dependencies on the two prey species, whereas the 
globalist have low values of them. 
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Figure 3: Habitat description of Example 5.0.1 



Hence if we simulate the localist and globalist as predator species 1 and 2 in our model, then we would 
assume s% > S2 and p± > P2- We may at first leave the other parameters the same for the two predator 
species. In this situation, the stability diagram is as in the Figure 3. 




(d) (e) (f) 



Figure 4: Stability diagrams for localist vs. globalist. The six cases have the common values of intrinsic rate of increse r\ = 
r 2 = 1.5 and carrying capacities K\ = K 2 = 10 5 , and the first four has the same conversion ratios 511 = (712 = 921 = 922 = 0.01. 
The other parameters for the first four are (a) pi = 0.06, p 2 = 0.01, s l = 0.09, s 2 = 0.01, (b) pi = 0.06, p 2 = 0.01, s 1 = 0.09, 
s 2 = 0.012, (c) pi = 0.03, p 2 = 0.01, si = 0.09, s 2 = 0.01, (d) pi = 0.06, p 2 = 0.01, si = 0.03, s 2 = 0.01. (e) and (f) have 
the same value of pi = 0.06, P2 = 0.01, s± = 0.09, s 2 = 0.01 but they differ on the conversion ratios : (e) gn = a 2 ± = 0.13, 
912 = <?22 = 0.01, (f) c?n = 921 = 0.01, 912 = 922 = 0.011. 

Recall the red dots plotted on the diagram shows the tuples (an, 022) for which the system has positive 
and stable equilibrium. Hence if the system goes on stable state of coexistence, the tuple (011,022) must 
be close to the red region. As one can see, the range of 022 is very limited around 0.5, which means, the 
globalist must be generalist, that is, it should not be selective on the two local preys. On the other hand, 
the values of an are apart from the center 0.5. Hence the localist must be specialist, which means, selective 
on the two preys. Also consider the case (b) that generalist have slightly higher search rate S2 = 0.012 from 
0.01 of (a). As one can see, in (b) the red regions ranges for an have dramatically reduced down around 
the intervals (0.3,0.4) and (0.6,0.7). This means in this case the localist is not allowed to be too selective 
on preys. 
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(c) is the stability diagram for lower dependency si = 0.03 of localist. Now the localist must be more 
selective on preys, but the globalist can be little selective on preys. On the other hand, when the localist 
have lower search rate relative to its dependency, we have diagram (d), which shows pretty much restricted 
area of stable tuples (an, 0,22)- Thus in any cases, we conclude that the localist are forced to be specialist, 
and globalist to generalist, in order for long coexistence of all the four species. 

Now we consider other factors as well. First, suppose the localists are smaller than the globalists in body 
size. This difference will reflect in higher conversion ratios qn = (721 = 0.013 then the other (721 = 922 = 0.01. 
But one might guess that the diagram in this case would be similar to that for Figure 14, since increasing 
both qn and (721 woud have the same effect of increasing the search rate Si. See (e). Similary, the case 
of large localist described by increased value of qi2 = 922 = 0.011, (f) which will have the same effect of 
increasing the search rate S2, will likely to show a similar diagram as Figure 15. Hence regarding the body 
sizes, the interesting case will be asymmetric conversion ratios, say qn = 0.01, q\2 — 0.011, (721 = 0.011 and 
922 = 0.01. 

Adapted preys. Assume that prey 1 is adapted to predator 1 and prey 2 to predator 2. We represent 
this situation by setting the adaptation coefficients an = a 2 2 = 0.4, a 12 = a 2 i = 1.5. We give four examples 
corresponding to different set of search rates s and predator dependency p. 




en eu 
(c) (d) 

Figure 5: Common parameters : r\ = T2 = 1.5, K\ = if 2 = 10 5 , qn = qi2 = 921 = 922 = 0.01, an = a22 = 0.4, 
ai2 = asi = 1.5. (a) si = s 2 = Pi = P2 = 0.02 (b) Sl = s 2 = p 2 = 0.02, Pl = 0.03 (c) si = 0.07, s 2 = pi = p 2 = 0.02 (d) 
si = 0.07, pi = 0.04, s 2 = p 2 = 0.02. 

In (a), all the parameters other than the adaptation coefficients are the symmetric in species 1 and 2. 
Hence the diagram is symmetric about the diagonal y = x. Ignoring the narrow strips near the line en = 0.8 
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and e22 = 0.8, the dominating upper right region shows that in order for long persistence of the system, 
predator i should be far more efficient to prey i, which is much hard to catch than the other. Now we 
increase the dependency of predator 1 by 0.01 in (b). Now the vertical narrow strip is thicker, and the 
horizontal strip has been disappeared. As en < 0.8, predator 2 has wide variety of options for its hunting 
strategy but not allowed to be too effective on hunting prey 2. On the other hand, as en > 0.8, which 
means predator 1 is pretty much specialized for prey 1, then predator 2 must also be specialized for prey 2 
in similar extent. 

Instead, if we increase the search rate s\ from 0.02 to 0.07 in example (a), we get example (c). The 
upper-right region is the same, but we have upper-left region instead of the strap near the line en = 0.8 in 
example (b). Hence a long persistence will be also possible if both the predator 1 and 2 are specialized for 
prey 2. 

Lastly, example (d) shows that increasing the value of p\ in example (c) enlarges the upper left region, 
allowing wide range of options for predator 1. Comparing this example with (b), one sees that even if the 
predator 1 is more dependent on species, the coexistence condition might be more strict to the predator 2 
if predator 1 has high search rate. 

Small predators must be selective. Consider an ecosystem where two predators with different body 
size and two preys are competing. We want to investigate the stable tuples (an, 022) for which the system can 
persist with all the four species coexisting. To see this, set gn = 021 = 0.04 and 021 = g22 = 0.01, reflecting 
the predator specie 1 has small body size so that it has higher conversion ratio. Set other parameters 
symmetric for species 1 and 2. The stability diagram for this case is dipicted in Figure 5- (a), which looks 
like Figure 4-(a). As one can see, the stable values for an are apart from the neutral 0.5; hence the small 
predator, specie 1, must be selective on preys. On the other hand, the other predator can be neutral on 
preys, but they are forced to prefer the other preys than those the small predator prefers. In fact, as 
we pointed out in previous examples, high gn,g2i values does the same for high s\ values. Hence small 
predators would be similar to the one who searchs preys more. 

More dependent on preys, more selective on preys. Set p\ = 0.05 andp 2 = 0.01, so that predator 
1 is five times more dependent on the preys. The stability diagram Figure 5-(b) shows that Predator 1 must 
be selective on preys to in order for coexsistence. 

Assymetric diagrams. The stability diagrams we have observed so far were all symetric about the 
point (0.5,0.5), which means that there is no difference between the roll of species 1 and 2. Now let us 
consider more complex combination of parameters which yields assymetric stability diagrams. The diagram 
Figure 5-(c) shows predator 2 must be either specialized for prey specie 2 or neutral, while predator 1 
can have any value of an, depending on a22-The key factor that makes the diagram asymmetric in above 
example is the asymmetric values of conversion ratios g,j. Such combinations of conversion ratios indicates 
that there are some non-trivial interrelationship between the four species other than just the body sizes. If 
the conversion ratio were dependent only on body sizes, for instance, we can assume prey 1 have body size 1 
unit, so that predator 1 has size 100 unit since gn = 0.01. Similarly we would have that predator 2 has size 
50 = 1/(721 unit, and hence prey 2 has size unit 3 = 100 • a 2 i. However, the condition g 2 2 = 0.01 indicates 
that predator 2 has hundred times of body size of that of prey 2, which is not the case here. This is because 
the system of equations z\ = az 2l z 2 = bz 3 , z 3 = cz 4 and Z4 = dz\ may not have a solution in general. 

Figure 5-(d) is another example of such assymetric conversion ratios, with other parameters varied. This 
example describes a system where predator 1 searches more than predator 2, while the former are more 
dependent on preys, and prey 2 are more reproductive, and conversion ratios are assymetric. As depicted 
in the diagram, the stable tuples are pretty much restricted than early examples. 

Lastly, Figure 5-(e) shows that such combination of parameters actually determines the stretagics of 
predators. 
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(C) (d) (e) 



Figure 6: All the five cases have the same carrying capacities K\ = K% = 10 5 , and the first three has the same value of 
ri = T2 = 1.5 and s± = S2 = 0.01. The other parameters are given as (a) pi = P2 = 0.01, gn = q2i = 0.04, q±2 = 922 = 0.01, 
(b) pi = 0.05, p 2 = 0.01, (jn = gi2 = 921 = 922 = 0.01, (c) pi = 0.03, p 2 = 0.01, qn = 0.01, 912 = 0.03, 021 = 0.02, 922 = 0.01. 
(d) has ri = 1.2, r 2 = 1.5, pi = 0.03, p 2 = 0.01, 911 = 0.01,912 = 0.015,921 = 0.02, 922 = 0.01, si = 0.02, s 2 = 0.01 and (e) 
has ri = 1.3, r 2 = 1.7, pi = p 2 = 0.01, 911 = 0.03, 912=0.01 , 921 = 0.015, g 22 = 0.02, Sl = 0.04, s 2 = 0.01. 

5. Concluding remarks 

6. Appendix 

6.1. comparison of the four- dimensional discrete-time Lotka-Volterra system to the corresponding continuous- 
time system 

6.2. comparison of the four- dimensional discrete-time Lotka-Volterra system to the discrete- system without 
self-inhibition term for prey species 

6.3. definitions of oscilators and basin of oscillation 
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